#!/usr/bin/perl -w
#-----------------------------------------------------------+
#                                                           |
# Fasta2Dir.pl                                              |
#                                                           |
#-----------------------------------------------------------+
#                                                           |
#  AUTHOR: James C. Estill                                  |
# CONTACT: JamesEstill_@_gmail.com                          |
# STARTED: 08/15/2006                                       |
# UPDATED: 08/15/2006                                       |
#                                                           |
# DESCRIPTION:                                              |
#  For a given fasta file, generates a file that describes  |
#  some information about the sequences including length.   |
#  This will be tab delim so that the length data can be    |
#  used to created a histogram in R. The main use for this  |
#  program will be to create a histogram of gene lengths    |
#  and TE lengths for Rice and Arabidopsis to compare to    |
#  the average sequence read length to illustrate why       |
#  assembly is necessary.                                   |
#                                                           |
#-----------------------------------------------------------+

print "The Get Seq Size program has started\n";
#exit;

use Bio::SeqIO;
use Cwd;

# Get command-line arguments, or die with a usage statement
my $thresh = '700'; # Threshold length
my $tcount = '0';   # Number of fragments bigger then threshold
my $usage = "Fasta2Dir.pl infile outfile \n";
my $infile = shift or die $usage;
my $outfile = shift or die $usage;
my $CurrentDir = cwd();                      
my $OutputDir = $CurrentDir."/CountOut/";
my $OutPath = $OutputDir.$outfile;
mkdir $OutputDir;

$SeqNum = 1;

# create one SeqIO object to read in,and another to write out
my $seq_in = Bio::SeqIO->new('-file' => "<$infile",
                             '-format' => 'fasta' ) ||
    die "Could not open the file ".$infile."\n";

# Open outfile to write output to
open (OUT, ">".$OutPath) ||
    die "Could not open output file:".$OutPath."\n"; 

# write each entry in the input file to the new output file
while (my $inseq = $seq_in->next_seq)
{
    $SeqNum++;


    my $SeqId = $inseq->primary_id();
    my $SeqLen = $inseq->length();

    # INCREMENT THRESHOLD IF NECESSARY
    if ($SeqLen > $thresh){ $tcount++; }

    # PRINT THE INFO TO THE OUTPUT FILE
    print OUT $SeqId."\t".$SeqLen."\n";

    # PRINT TO SCREEN TO SHOW WERE WE ARE
    print $SeqNum."|".$SeqLen."\n";;
}

# CLOSE THE OUTPUT FILE
close OUT;

#-----------------------------+
# PRINT SUMMAR TO THE SCREEN  |
#-----------------------------+
print $SeqNum." sequences have been processed\n";
print $tcount." > ".$thresh."\n";
my $perover = ($tcount/$SeqNum)*100;
print "percent over thres is :".$perover."\n";

exit;
